// test_eign_calc.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include <iostream>
#include <Eigen/Dense>
//using Eigen::MatrixXd;
using namespace Eigen;
using namespace std;
void test_addition_subtraction()
{
	Matrix2d a;
	a << 1, 2, 
		 3, 4;

	MatrixXd b(2, 2);
	b << 2, 3,
		 1, 4;
	cout << "a = \n" << a << endl;
	cout << "b = \n" << b << endl;
	cout << "a + b = \n" << a + b << endl;
	cout << "a - b = \n" << a - b << endl;
	cout << "Doing a += b\n";
	a += b;

	Vector3d v(1, 2, 3);
	Vector3d w(1, 0, 0);
	cout << "-v + w - v = \n" << -v + w - v << endl;
	cout << "Now a = \n" << a << endl;
}
void test_scalar_multiplication_and_division()
{
	Matrix2d a;
	a << 1, 2,
		3, 4;
	cout << "a = \n" << a << endl << endl;

	Vector3d v(1, 2, 3);
	cout << "v = \n" << v << endl << endl;
	cout << "a * 2.5 = \n" << a * 2.5 << endl << endl;

	cout << "0.1 * v = \n" << 0.1 * v << endl << endl;

	cout << "doing v*= 2" << endl;
	v *= 2;
	cout << "Now v = \n" << v << endl << endl;
}

void test_transposition_and_conjugation()
{
	Matrix2d a;
	a << 1, 2,
		3, 4;
	cout << "a = \n" << a << endl << endl;

	cout << "a^T = \n" << a.transpose() << endl << endl;

	cout << "a conjugate = \n" << a.conjugate() << endl << endl;
	cout << "a * = \n" << a.adjoint() << endl << endl;
}

void test_matrix_matrix_and_matrix_vector_multiplication()
{
	MatrixXd mat(3,2);
	mat << 1, 2,
		   3, 4,
		   5, 6;

	Vector2d u(-1, 1), v(2, 0);
	cout << "mat * mat = \n" << mat.transpose() * mat << endl << endl;
	cout << "mat * u = \n" << mat * u << endl << endl;
}
void test_dot_product_and_cross_product()
{
	Vector3d v(1, 2, 0);
	Vector3d w(0, 1, 0);

	cout << "v.dot(w) = \n" << v.dot(w) << endl << endl;

	cout << "Cross product:\n" << v.cross(w) << endl << endl;
}

void test_basic_arithmetic_reduction_operations()
{
	Matrix2d mat;
	mat << 1, 2,
		3, 4;

	cout << "mat.sum = " << mat.sum() << endl << endl;
	cout << "mat.prod = " << mat.prod() << endl << endl;
	cout << "mat.minCoeff = " << mat.minCoeff() << endl << endl;
	cout << "mat.maxCoeff = " << mat.maxCoeff() << endl << endl;
	cout << "mat.trace  = " << mat.trace() << endl << endl;
}

void test_block_operator()
{
	Matrix4f m;
	m << 1, 2, 3, 4,
		5, 6, 7, 8,
		9, 10, 11, 12,
		13, 14, 15, 16;

	cout << " m = \n" << m << endl << endl;

	//m.topLeftCorner(1, 3) = m.bottomRightCorner(3, 1).transpose();
	m.block<1, 3>(0, 0) = m.block<3, 1>(m.rows()- 3, m.cols()-1).transpose();
	cout << "m = \n" << m << endl << endl;
}

void test_partial_reductions()
{
	MatrixXf m(2, 4);
	m << 1, 2, 6, 9,
		3, 1, 7, 2;
	MatrixXf::Index maxIndex;
	float maxNorm = m.colwise().sum().maxCoeff(&maxIndex);

	std::cout << "m = \n" << m << std::endl;

	std::cout << "col sum max Index: " << maxIndex << std::endl;
	std::cout << "col sum max Col : \n" << m.col(maxIndex) << std::endl;
	std::cout << "col sum = " << maxNorm << std::endl;
}
void test_linear_solving()
{
	Matrix3d A;
	A << 1, 2, 3, 4, 5, 6, 7, 8, 10;
	Vector3d b;
	b << 3, 3, 4;

	std::cout << "A = " << A << std::endl;
	Vector3d x = A.colPivHouseholderQr().solve(b);
	std::cout << "x = " << x << std::endl;

	Matrix2d B, b2;
	B << 2, -1, -1, 3;
	b2 << 1, 2, 3, 1;
	std::cout << "B = " << B << std::endl;
	Matrix2d x2 = B.ldlt().solve(b2);
	std::cout << "x = " << x2 << std::endl;
}

void test_relative_error()
{
	MatrixXd A = MatrixXd::Random(100, 100);
	MatrixXd b = MatrixXd::Random(100, 50);
	MatrixXd x = A.fullPivLu().solve(b);
	double relative_error = (A * x - b).norm() / b.norm();
	std::cout << "The relative error is : \n" << relative_error << std::endl;
}
void test_eigen_value_and_eigen_vector()
{
	Matrix2d A;
	A << 1, 2, 2, 3;
	std::cout << "A = \n" << A << std::endl;

	SelfAdjointEigenSolver<Matrix2d> eigenSolver(A);
	if (eigenSolver.info() != Success)
	{
		std::cout << "calc eigen values failed!\n";
	}

	std::cout << "eigen values = \n" << eigenSolver.eigenvalues() << std::endl;
	std::cout << "eigen vectors = \n" << eigenSolver.eigenvectors() << std::endl;
}
int main()
{
	/*MatrixXd m(2, 2);
	m(0, 0) = 3;
	m(1, 0) = 2.5;
	m(0, 1) = -1;
	m(1, 1) = m(1, 0) + m(0, 1);
	std::cout << m << std::endl;*/

	/*Matrix3d m;
	m << 1, 2, 3,4, 5, 6,7, 8, 9;
	std::cout << std::endl << "m1 = " << m << std::endl;
	m += Matrix3d::Constant(3, 3, 1.2) * 50;

	std::cout << "m2 = " << m << std::endl;

	Vector3d v(3);
	v << 1, 2, 3;
	std::cout << "m * v = " << std::endl << m * v << std::endl;*/
	//test_addition_subtraction();
	//test_scalar_multiplication_and_division();
	//test_transposition_and_conjugation();
	//test_matrix_matrix_and_matrix_vector_multiplication();
	//test_dot_product_and_cross_product();
	//test_basic_arithmetic_reduction_operations();
	//test_block_operator();
	//test_partial_reductions();
	//test_linear_solving();
	//test_relative_error();
	test_eigen_value_and_eigen_vector();
	system("pause");
    return 0;
}

